Mon. Not. R. Astron. Soc. 000, QHH] (2013) Printed 10 April 2013 (MN LMftX style file v2.2) 



The Warm DM halo mass function below the cut-off scale. 



Raul E. Angulo 1 *, Oliver Hahn 2 f, Tom Abel 1 :):. 

1 Kavli Institute for Particle Astrophysics and Cosmology, 

Stanford University, SLAC National Accelerator Laboratory, Menlo Park, CA 94025, USA 

2 Department of Physics, ETH Zurich, CH-8093 Zurich, Switzerland. 



CO 

o 

<N 

Oh' 
< 
00 



o 

U 

Oh 

6 

H— » 

a 



> 

o 

o 
m 



10 April 2013 



% 



ABSTRACT 

Warm Dark Matter (WDM) cosmologies are a viable alternative to the Cold Dark Mat- 
ter (CDM) scenario. Unfortunately, an accurate scrutiny of the WDM predictions with 
N-body simulations has proven difficult due to numerical artefacts. Here, we report on 
cosmological simulations that, for the first time, are devoid of those problems, and thus, 
are able to accurately resolve the WDM halo mass function well below the cut-off. We 
discover a complex picture, with perturbations at different evolutionary stages populating 
different ranges in the halo mass function. On the smallest mass scales we can resolve, 
identified objects are typically centres of filaments that are starting to collapse. On inter- 
mediate mass scales, objects typically correspond to fluctuations that have collapsed and 
are in the process of relaxation, whereas the high mass end is dominated by objects similar 
to haloes identified in CDM simulations. When explicitly show how the formation of low- 
mass haloes is suppressed, which translates into a strong cut-off in the halo mass function. 
This disfavours some analytic formulations that predict a halo mass function that would 
extend well below the free streaming mass. We argue for a more detailed exploration of 
the formation of the smallest structures expected to form in a given cosmology, which, we 
foresee, will advance our overall understanding of structure formation. 

Key words: cosmology: theory - large-scale structure of Universe. 



1 INTRODUCTION 

Since the 1980s, Warm Dark Matter (WDM) has been an at- 
tractive alternative to Cold Dark Matter (CDM) as the main 
gravitating component in the Universe. For a long time, how- 
ever, WDM was disfavoured compared to CDM due to, in 
part, the additional free parameter required (the mass of the 
DM particle). This has changed recently, and WDM has again 
attracted the attention of the cosmological community as a 
viable and competitive cosmological model. Traditionally, the 
gravitino was favou red as the hypothe tical particle that could 
serve as WDM fe.g. lMoroi et al.l ll993, and references therein), 
while, more rece ntly, the interest has f ocused on the sterile 
neutrino (see e.g. iBovarskv et al.ll2009bl . for a recent review). 
The reason for the recent revitalisation of WDM is that 
a suitable AWDM model could reproduce all the successes 
of ACDM on large scales, and, in addition, it may alleviate 
the alleged tension between the CDM model and some ob- 
servations. Such measurements concern mainly mass scales of 
~ 10 10 Mq and include: the dynamics of Milky- Way satellites 
(|Bovlan-Kolchin et al.ll201ll : lLovell et alll2012h . the velocity 



function of HI selected galaxies (jZavala et al.ll2009h . and the 
abundance of low mass galaxies at low and high redshifts 
(|Menci et alJl2012h . 

The key feature of WDM cosmologies that distinguishes 
them from CDM, is the lack of initial small-scale density 
fluctuations. Before recombination, WDM particles have rel- 
atively high velocities (set when these particles become non- 
relativistic) , which allows them to travel further than a "free- 
streaming" distance Rf s . Thus, particles move out of over- 
dense regions of size Rf s and smaller, and therefore inflation- 
generated small-scale density and potential perturbations are 
washed out. This dissipation is captured as a strong sup- 
pre ssion of the mass tran sfer function below a "cut-off" scale 
(cf. lBond fe Szalavlll983h . This difference relative to CDM is 
expected to affect the abundance of collaps ed objects (DM 



haloes) below the cut- off mass scale 
120121 : iMenci et al.ll2012l : iBenson et al 
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20 



ISchneider et al.l 
and also the in- 



ternal properties of haloes (| Colin et al.l 1 2000, e.g.). In par- 
ticular, it is expected the abundance of low-mass haloes to 
be suppressed, and the halo density profiles to be less con- 
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1 Note that the free-streaming scale today is much smaller than that 
at early times. Thus, scales below the cut-off scale in the transfer 
function are, in principle, able to collapse gravitationally, but not 
those below today's value of the free-streaming scale. 
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centrated. In addition, the formation time, as well as the ac- 
tual halo formation mechanism is modified, as a larger frac- 
tion of haloes form directly from the collapse of filaments. 
Additionally, phase correlations of the overall density field 
are expected to be rela ted to the mass of the DM particle 
(lObreschkow et al.ll20l3 ). 

Besides the obvious cosmological interest in WDM sce- 
narios, this class of models also offers an interesting test case 
from a theoretical point of view, especially for cosmological 
simulations. Free- streaming provides a small-scale limit to the 
structure formation problem, thus, in principle, it is possible 
to capture the full hierarchy of objects expected to be relevant 
during the formation of a halo. This can provide clues as to 
how hierarchical structure formation proceeds in general, and 
its connection to cosmological parameters and the spectrum 
of density fluctuations. In contrast, this is less direct in CDM 
since there are always structures and interactions not resolved 
by any given mass and force resolution at late times. 

Unfortunately, an accurate numerical scrutiny of the pre- 
dictions of a WDM model has proven difficult. Despite a mul- 
titude of attempts, N-body simulations, and analytical frame- 
works, have not been able to decisively quantify even one the 
simplest properties of the nonlinear field: the abundance of 
DM haloes. Cosmological WDM simulations show the dom- 
inant presence of an artificial population of low-mass DM 
haloes. This phenomenon originates in the numerical (unphysi- 
cal) fragmentation of filaments, and it exceeds by far the popu- 
lati on of real WDM haloes on scales near the cut-off scale (see 
eg. 



Avila-Reese et al.ll200ll ; [Bode et alJl200lt IWang fc White] 



120071 ; lMelottll2Q07h . Although high er mass resolution hel ps to 
somewhat reduce this problem fe.g. lSchneider et al.l2013n . the 
mass scale dominated by artificial fragmentation changes only 
weakly with an increase of the mass resolution. Besides the 
obvious impact of artificial fragmentation in the halo mass 
function, it is possible that they also affect the internal prop- 
erties of DM haloes, which grow by accreting these objects. 
Finally, it is currently unknown what type of structures, if 
any, inhabit below the cut-off mass scale. The latter is rele- 
vant when modelling galaxy formation in WDM scenarios. 

Similarly, recent analytic approaches, based on 
modifications of the Extended Pre ss- Schechter theory 
(jPress fe Schechterl H97I lLacev fe Colel Il993h . produce pre- 
dictions that differ by orders of magnitude from each other 
in the way in which t he halo mass function is suppressed in 
WDM (compare e.g. ISmith fc Markovid I20T1I : iBenson et all 
I2013T ). The reason for this is mainly uncertainties in the 
formulation of traditional excursion sets, together with 
possible modifications to the shape of the barrier for collapse 
near the cut-off scale. 

In this paper, we revisit the issue of the halo mass func- 
tion in WDM cosmologies and propose an answer to the 
questions posed above. In order to overcome previous lim- 
itations, we carry out a suite of high-resolution cosmolog- 
ical N-body simulations that feature a recently develope d 
method to compute gravitational forces ([Hahn et al.l 120120 . 
This new technique allows us to alleviate some of the long- 
known problems originating from employing excessive force 
resolutio n compared to the mass re s olution in the simu- 
lations (lEfstathiou h Eastwoodl Il98ll: ICentrella et aD Il988l : 
iMelott fc Shandarinlll989l : ISplinter et al.l Il998h . and also to 
reduce discreteness noise in the large-scale density and tidal 
fields, which we find are key players in causing artificial frag- 
mentation. 



With these tools at hand, we are able to robustly compute 
for the first time the WDM halo mass function at, and be- 
low, the cut-off mass scale. Our simulations unveil systems of 
different characteristics, and at different stages of formation, 
populating different mass ranges in the halo mass function. 
Explicitly, well below the cut-off scale we find dense filaments 
and sheets that have started to collapse into 3D systems. At 
larger mass scales we find proto- haloes, systems that are col- 
lapsed but still in the process of virialisation. Only above the 
cut-off scale we observe systems traditionally regarded as DM 
haloes. When we consider only halos, We observe a strong sup- 
pression of the halo mass function, together with a cut-off on 
small masses, however. However, due to the wealth of different 
types of structures, the position of the cut-off depends on the 
actual halo definition one wishes to adopt. 

The structure of the paper is as follows: In §2 we provide 
details of the simulation techniques and the construction of the 
halo catalogues. In §3 we present our results, focusing on the 
abundance of collapsed objects, the WDM halo mass function, 
and exploring the characteristics of objects located at different 
mass ranges in the halo mass function. We then discuss our 
findings and possible implications. Finally, in §4 we provide a 
summary of our work along with directions for future work. 



2 NUMERICAL TOOLS 

In this section, we describe the numerical simulations we use to 
study dark matter haloes in a WDM scenario. We also describe 
our halo identification procedure and the construction of the 
group catalogues. 



2.1 Initial conditions 

We start by computing the power s pectrum of densi t y fluc - 
tuations using the fitting formulae of lEisenstein fc Hul ([19991 ). 



We adopt a set of cosmological parameters consistent with 
the published measu rements of the WMAP7 data release 
(iKomatsuet al.ll2010h . Explicitly: Q m = 0.276, Ov = 0.724, 
n b = 0.045, h = 0.703, a 8 = 0.811 and spectral index 
n s — 0.96. Note we set the normalisation of the power spec- 
trum according to as via a CDM spectrum, so that the am- 
plitude of fluctuations on large-scales is independent of the 
WDM particle mass. 

We then incorporate the effects of a thermally produced 
warm dark matter pa rticle on the trans fer function following 
the fitting formula of lBode et all (J200lh Fl Explicitly, a fit to 
the WDM transfer function of density perturbations is given 
by: 



T W DM(k) = T CDM (fc) [1 + (ak) 2 



(1) 



with 



2 See also lViel et al.l (120051 ) who adopt a very similar fitting func- 
tion but indicate a simple way in which (non-thermal) sterile neutri- 
nos can be accounted for as well by a change to an effective thermal 
mass (cf. also lColombi et aLl llQQG). 
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where radm is the DM particle mass (or the effective sterile 
neutrino mass), in units of keV, and gx is the number of de- 
grees of freedom that the particle contributes to the number 
density which is 3/2 in our case. 

In this paper we will explore the case where mam = 
250 eV. Our choice is inconsistent with current constraints 
placed by observations of the Ly-a forest power spectrum 
(which set a lower lim it at the order of keV ()Viel et al.l 12005k 
iBovarskv et al.ll2009al )). However, a low md m has the advan- 
tage of allowing us to resolve structures at redshift zero much 
below the cutoff mass scale using only relatively modest mass 
resolution and computational resources. In addition, our re- 
sults can be readily extended and generalised to others DM 
particle masses. 

For our choice of cosmological parameters and WDM par- 
ticle mass: a = 0.26/i _1 Mpc, which translates into a free- 
streaming mass-scale 



M ft 



47T 



p(a/2) 3 ~7x 10 8 /i _1 M G 



(3) 



and a half-mode mass-scale (cf. I Colin et al.l 120081 ; 
ISchneider et~aTll2012h 



M hm ~ 4.3 x 10 3 Mts ^ 3.0 x 10 12 /i _1 M . 



(4) 



This mass-scale is where a deficit in the mass-function com- 
pared to CDM is expected to be a factor of two. 

Using the WDM pri mordial power spec trum discussed 
above and the MUSIC code ([Hahn &; AbefeOllh , we then create 
the initial particle configuration at z — 63 for our numerical 
experiments. We do this by perturbing a particle distribution, 
initially arranged in a regular cubic lattice, according to the 
Zel'dovich approximation. 

We note that we have not attempted to explicitly include 
thermal velocities (on top of gravitationally induced veloci- 
ties) in our initial conditions, since it is both negligible for 
the results of this paper as well as unclear how a proper im- 
plementation would proceed. The RMS velocity of the micro- 
scopic WDM particles following the Fermi-Dirac distribution 
(for t hermally produced warm dark matter) is given at redshift 
z by (iBode et al]l200lh 



0.043 (1 + z) 

1/3 



(5) 

fej fej UJ fej kms • 

For the radm = 0.25 keV chosen in this paper, we find 
a v — 0.28 kms -1 at z = 0. Clearly, a value considerably 
smaller than those arising due to nonlinear structure forma- 
tion. For instance, the virial velocity of a typical DM halo 
in our simulations is 100-1000 kms -1 . Of course, this RMS 
velocity can be of relevance in the most central parts of the 
halo es, determining details of the phase-space density there 
(e.g. iDalcanton fc H ogan 20 011V as well as the thickness of 
caustics (e.g. IWhite &; Vogelsbergerll2009h . 

Although the RMS velocity corresponds to a microscopic 
value, it is sometimes regarded as a macroscopic one, and im- 
plemented in N-bo dy simulations as random kicks of sim ula- 
tion particles (e.g. IColm et al.ll2008l : iMaccio et alJl2013h . We 



emphasise, however, that this is simply an Ansatz and that 
simulation particles represent a coarse-grained phase-space 
distribution, thus each of them averages over a statistical en- 
semble with a negligibly small dispersion around this mean. 
On the contrary, a kick to a single simulation particle is equiva- 
lent to a locally coherent motion of a large ensemble of actual 
WDM particles. This leads to a velocity spectr um inconsis- 
tent w ith the results of linear perturbation theory (| Colin et al.l 
I2008T ) , and it is clear that the evolution of this numerical setup 
is not equivalent to the evolution of the fine-grained distribu- 
tion function. Therefore, neglecting this dispersion is a very 
good, as well as convenient, approximation when one is only 
concerned with the mass of DM haloes, as it is in our case. 



2.2 Gravitational Evolution 

We perform a series of cosmological N-body simulations evolv- 
ing 1024 3 particles inside a cubical region of L = 80 h~ 1 Mpc a 
side. For our choice of cosmological parameters, each of these 
simulation particles has a mass equal to 3.65 x 1O 7 /i -1 M . 
This mass resolution and volume is sufficient to have a fair 
sample of haloes located at the half-mode mass, which is re- 
solved with almost 100, 000 particles. 

We evolve simulation particl es using a mem ory- efficient 
version of the P-Gadget3 code ([Springell I2005T L which was 
originall y developed and op timised for the Millennium-XXL 
project (|Angulo et al.l I2012T ). In this code we have imple- 
mented three different methods to compute gravitational 
forces. In the remainder of the paper we refer to them as 
Tree-PM, PM and T4PM, and are described in the following: 

1. Tree-PM: This method corresponds to the standard nu- 
merical configuration followed in state-of-the-art calculations. 
Long-range interactions are calculated using a PM method 
(|Hocknev fc Eastwood! 1 1 98 lh . whereas short-range forces are 
calculated using a multipole e xpansion of the forc e field to- 
gether with a Tree algorithm ([Barnes & Hutl ll986). In order 
to reduce two-body scattering and binary particle systems 
(among other artefacts), forces need to be softened on small 
scales. We do this in our runs using a Plummer-equivalent 
softening length equal to 5/i -1 kpc. 

2. PM: Here, gravitational forces are only given by the PM 
method, i.e. we compute the gravitational potential field on a 
grid by solving the Poisson equation using Fourier methods. 
In our runs we use a grid of 2048 3 points and forces are Gaus- 
sian smoothed on scales roughly equal to twice the grid size, 
2 Ax = 80/i -1 kpc. This length scal e matches the mean inter- 
particle separation. As argued by lAngulo et al.l (|2013r ), this 
numerical configuration suppresses undesired collisionality of 
the N-body system, and it is particularly successful (compared 
to a Tree-PM run) in following accurately the gravitational in- 
teraction of baryons and DM. 

3. T4PM: T his method is an im plementation of the algorithm 
proposed bv lHahn et al.l (|2012h . In short, a Delaunay triangu- 
lation of the Lagrangian particle distribution defines a phase- 
space element (a tetrahedron) that can be reconstructed at 
any desired later time to reconstruct the respective density 
field. At all times, the density of each tetrahedron is assumed 
to be uniform. In practice, we represent the contribution of 
each tetrahedron to the total mass field using 4 virtual par- 
ticles (they carry mass, but do not interact directly with the 
fluid) whose spatial distribution matches the monopole and 
the quadrupole of the parent tetrahedra. We deposit these 
particles onto a 4096 3 mesh with CIC interpolation and com- 
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Figure 1. Images of the cosmic density field at z = in a WDM 
cosmology with m^m = 0.25/ceV, as predicted by the three dif- 
ferent methods to compute gravitational forces. Each image dis- 
plays a projection of a 20 h~ 1 Mpc thick slab, of a region of size 
80 x 40 /i _1 Mpc, which focuses on the most massive halo present in 
the simulation volume. Additionally, in the three panels all haloes 
with mass M200 > 1 X 10 10 /i _1 M0 are displayed using red circles 
whose radii equal the haloes' virial radii, i?200- Note that the dif- 
ferent methods show different amount of spurious fragmentation, 
which results in sequences of haloes aligned along filaments. 



pute forces using the PM method, as described above. The 
spatial resolution of these runs is 40/i _1 kpc, twice as high as 
in the PM case. 



In Fig. [T] we show the density projection of a subregion 
of our WDM simulations, as computed by the three methods 
described above. We will analyse this figure in detail in the 
next section, but it is readily apparent that all three runs 
display the same large-scale structure, but they differ on small 
scales. 

Additionally, we have carried out a set of CDM N-body 
simulations, with identical force and mass resolution as the 
WDM runs. The initial conditions of all simulation boxes were 
generated using identical white-noise fields, which simplifies 
their comparison. 



2.2.1 Computational Performance 

The computational resources required for the three methods 
are similar but differ systematically. For our WDM runs, the 
Tree-PM, PM and T4PM runs require about 3000, 1000 and 3500 
CPU-hours, respectively. For the Tree-PM run, ~ 50% of the 
time was employed in the construction and walking of the 
tree, while the overhead of the T4PM respect to the PM method 
is caused by our on-the-fly calculation of the initial tessella- 
tion together with the position calculation and depositing of 
the mass-carriers particles. The peak memory consumption 
of the T4PM run is ~ 600Gb, to be compared with ~ 150Gb 
employed by the PM run. The difference is dominated by an 
extra set of pointers needed to recover the initial connectiv- 
ity of the phase-space tessellation. We remark that the extra 
factor of ~ 4 is small considering that we effectively have 24 
times more particleqj representing the density and force field. 
For a dramatic increase in force accuracy it hence seems very 
worthwhile to afford this additional cost in memory and run 
time. In addition, our implementation is suboptimal in terms 
of memory consumption: it is possible to carry out the T4PM 
run with a memory footprint identical to that of PM, at the 
cost of slightly more CPU time. 



2.3 Dark Matter Haloes 

For each of our simulations, w e produce on-the-fly friends-of- 
friends fFoF. lDavis et al.lll985T ) halo catalogues. We use a non- 
standard linking length parameter of b = 0.05 times the mean 
inter-particle separation, keeping objects with 20 or more par- 
ticles. This unusual choice of b (as compared to b = 0.2) is 
required to avoid large FoF haloes percolating the cosmic web. 
We will return to this point below. 

For each FoF halo, we compute a spherical-overdensity 
(SO) mass, taking the centre of mass of the parent FoF group 
as the SO centre. We define the halo boundary as the sphere 
of radius R200, which contains a mean density of 200 times 
the critical density, pcrit- Therefore, the mass of the halo is 

M20O = |7Ti^ioO 200 pern. 

We discard substructures from our catalogues whose R200 
spheres overlap with that of a more massive halo. At z — 
this procedure finds 8359, 3422 and 2916 objects with mass 
M200 > lO lo /i _1 M in the WDM Tree-PM, PM, and T4PM runs, 
respectively. These are a factor of 15 — 40 smaller than in a 
CDM Tree-PM run, where we detect 127' 133 structures. In 
Fig. [TJ we overplot this halo catalogue on top of the dark 
matter density field. 

In passing, we note that even though the FoF algorithm 
with the standard choice, b = 0.2, works satisfactorily for the 
Tree-PM run, it fails to deliver a reasonable halo catalogue for 
the other two cases. While in the Tree-PM case the filaments 
are broken into small haloes, in the T4PM and PM runs a strong 
artificial fragmentation is absent (as can be seen in Fig. [TJ, 
thus filaments and sheets are more homogeneous with sharp 
dense cores: genuine two and one dimensional dense structures 
exist. The FoF algorithm links these filaments to dark matter 
haloes located at their ends, and with other nearby haloes. We 
show examples of this problem in Fig. [2] where we display a 
projection of the particles associated with the two most mas- 
sive FoF haloes at z = 1 (top row) and z = (bottom row), 



3 Each particle contributes to six distinct tetrahedra, and each 
tetrahedra is represented by four particles. 
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Figure 2. The 'standard' linking length b = 0.2 selects large 
parts of a WDM simulation once the forces are captured accu- 
rately enough that filaments do not artificially fragment. Den- 
sity projections of the particles belonging to the two most mas- 
sive FOF- "haloes" in our WDM T4PM simulation of a 250eV dark 
matter model. Objects at z = 1 (top row) and z = (bottom 
row) are shown. These haloes have a mass of 2.7 x 1O 14 /i -1 M0 
and 1.6 x 10 14 h^M® (z = 1), and of 6.4 x 10 14 /i- 1 M o and 
2.8 x lO^-iM© (z = 0). 



as identified in the T4PM runs. The largest b — 0.2 FoF struc- 
ture at z — has a mass of 6 x 10 14 h~ 1 MQ and spans almost 
10/i _1 Mpc. At z = 1, the failure of the FoF algorithm is even 
worse: the biggest FoF halo spans almost one quarter of the 
simulation box size! 

In order to avoid such problems, we employed a small 
linking length that ensures that only local high density peaks 
are selected as the starting point for our SO halo catalogues. 
We tested that the resulting SO halo mass function was insen- 
sitive to small changes in b about our preferred value of 0.05. 
However, for values approaching b = 0.2, the mass function 
agrees only with the cases with smaller linking lengths at the 
highest mass end. On any other mass scale, it shows a noto- 
rious deficit of structures. This is because a large fraction of 
small haloes are artificially linked to form a single larger FoF 
structure, and thus are not present in our list of SO candidate 
haloes. 



3 RESULTS 

The main goal of this paper is to quantify the abundance of 
haloes expected in WDM cosmologies, especially below the 
cut-off scale. An accurate account of this is important, firstly, 
to establish robustly the predictions of WDM which can then 
be tested against observational data, and secondly, to under- 
stand more generally the collapse and assembly of DM haloes 
in the presence of a resolved cut-off scale in the perturbation 
spectrum. This in turn can help to understand the formation 
and properties of micro-haloes expected for some CDM par- 
ticle candidates (e.g. the neutralino). In addition, these cos- 



3.1 Halo abundance 
method 



dependence on the numerical 



Previous numerical simulations have not been able to ex- 
plore the cut-off mass scale because it is dominated by a 
population of low-mass haloes aligned within filaments. This 
phenomenon has been reported i n numerical simulation for 



decades: including early works (e.g.lMelott &: Shandarinll 19891 ; 
lAvila-Reese et al.ll200ll : lBode et alJl200ll:lKnebe et al.l | 2003h . 



and also recent state-of-the art runs (|Wang &; White! I2007J ; 
lLovell et alll2012l : ISchneider et alJl2012h . 

Initially, it was not clear whether a real and physical frag- 
mentation of filaments could be in place, or if it has its origin in 
numerical inaccuracies. This has settled recently, and there is a 
consensus that these haloes are numerical artefacts. Evidence 
for this is that their spatial distribution is closely related to 
the initial unperturbed particle load, and that their abundance 
changes (albei t slowly) with mass re solution (in fact oc m p ). 
In particular, IWang fc White! (|2007n have analysed this prob- 
lem in detail and concluded that their presence is caused by 
non-zero small-scale fluctuations of the ID-projected densit 
field. This is related to warnings of lMelott & Shandarinl (|198' 
about using excessive force resolution compared to the mass 
resolution as it leads to fragmentation also in 2D cosmological 
simulations. This is indeed the regime in which state-of-the-art 
simulations are carried out: the typical force resolution is set 
to a value 10 — 100 times smaller than the mean interparticle 
separation. 

The numerical nature of the fragmentation is also illus- 
trated in Fig. [TJ where we show a density projection of a 
20 /i -1 Mpc thick slab through our three simulation boxes at 
z = 0. The visualization technique is identical for all three 
panels, and corresponds to a CIC density (for the T4PM run, 
we project the flow tracers, not the 24 times more abun- 
dant mass carriers); thus any difference is a result of real 
discrepancies in the spatial distribution of particles. In these 
images, we overplot the halo distribution over the underly- 
ing dark matter field. We display only haloes with SO mass 
M200 > 2 x lO lo /i -1 M0, which is the resolution limit of our 
simulations, as we will discuss below. It is straightforward to 
see that all three runs, which use different methods to com- 
pute gravitational forces, display the same large-scale struc- 
ture while, however, differences exist on small scales. 

In the top panel, we display results obtained with the 
most commonly used method to compute gravitational forces 
(c.f.§2.1), labelled as Tree-PM. Fragmentation of filaments into 
small clumps is clearly visible in several places, for instance, in 
the two filaments located in the lower half of the image. These 
clumps are indeed very dense and are identified as haloes by 
our FoF- SO algorithm, and are thus highlighted by red cir- 
cles. In the middle panel, which shows the PM simulation, 
forces are effectively softened below the mean inter-particle 
separation and low-mass haloes aligned with the filaments are 
considerably less abundant. In the bottom image, displaying 
the T4PM run, artificial fragmentation virtually does not exist! 
Even though this run has a force resolution twice as high as 
in the PM case. 

Therefore, we see that the fragmentation of filaments is 
closely related to the force calculation, or more precisely, to 
the combination of force and mass resolution. We note that, 
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Figure 3. Comparison of the z = SO halo mass function in a 
WDM scenario for a 250eV dark matter particle mass. Coloured 
lines show the predictions of three different methods to compute 
gravitational forces: standard Tree-PM (red line), only PM (blue 
line), and the method of Hahn et al. (2012) (blue line) We also 
display the mass function expected for a CDM case for comparison 
(black line). Vertical dashed lines indicate a limit where the abun- 
dance of haloes is not affected by finite force resolution. The bottom 
panel shows the ratio of our results to the expectations for a CDM 
particle. 



in these simulations, the crucial difference is not the actual 
method to compute forces (e.g. a PM versus a Tree+PM), but 
the chosen force resolution for a given mass resolution. We 
have explicitly tested this assertion by varying the size of the 
mesh in the PM run. An excessive force resolution causes local 
minima in the global potential around simulation particles, 
which eventually grow and accrete neighbouring particles. In 
this sense the T4PM method has the advantage of smoothing 
these minima since it provides a smoothe r representation of 
the mass field (see e.g. iKaehler et al.ll2012i ) and thus more ac- 
curately captures the smooth but dense structure of the den- 
sity field in regions of strong anis otropic com pression. This was 
already qualitatively seen by lHahn et al.l (|2012h , who found 
that the T4PM method suppresses artificial fragmentation in 
WDM scenarios. The advantage of T4PM has the price that 
the density in the inner regions of haloes is overestimated. 
This is because the evolution of highly distorted Lagrangian 
phase-space elements in regions of strong mixing cannot be 
represented correctly by the piecewise linear approximation 
to the distribution function that is only tracked by the La- 
grangian motion of the particles. In principle this limitation 
can be overcome by an adaptive mass refinement (Hahn, An- 
gulo & Abel, 2013, in prep.). Nevertheless, this limitation has 
a very minor effect on the halo masses, and thus, on our results 
regarding the halo mass function. 

In Fig. [3l we can quantitatively see the differences in the 



predicted number of DM haloes at z — 0, as a function of 
their SO mass, M200, for our three methods. For comparison, 
we also display the halo mass function analogously constructed 
for a CDM simulation with matching volume and mass resolu- 
tion, and where forces are computed using the Tree-PM and PM 
method. For the latter, we use a mesh of 2048 3 cells, matching 
the force resolution of the T4PM run. The vertical dashed line 
indicates a mass of 2 x lO lo /i -1 M0, or, equivalently, ~ 700 
particles. This is an estimate for the mass limit above which 
we expect our results to be numerically robust. 

We choose this limit by comparing the resulting mass 
function in the CDM case for the PM and Tree-PM force meth- 
ods. Below M min = 2 x 10 10 /i _1 M©, the PM mass function 
shows a strong deficit of haloes, this is (1) because the force 
resolution of the PM run (80/i _1 kpc) is simply too low to re- 
solve (and keep bound) some low-mass density peaks, which in 
the Tree-PM run collapse as haloes; and (2) since it is plausible 
that even in CDM the lowest masses are a mixture of artifical 
fragments and true haloes, as suggested by the strong increase 
of the number of "peakle ss" haloes at low particle counts (cf. 
iLudlow fc Porcianill201lh . Above M m in this does not seem to 
be important and the mass functions are consistent with each 
other, showing only a small offset caused by a systematic un- 
derestimation of masses in the PM case. This is an indication 
that our results above M m [ n are numerically robust, and that 
the differences regarding artificial fragmentation are a result 
of our improved estimation of the force field and the respective 
reduction of discreteness effects, rather than due to a product 
of a somewhat low force resolution. Another aspect support- 
ing this is the fact that the amount of low-mass haloes in 
the T4PM run is lower than in the PM case, despite the former 
having higher force resolution. This is the opposite to what 
is expected if the suppression were caused by a lack of force 
resolution. 

As can be seen in the lower panel of Fig. [3l above M ~ 
3 x 10 11 h~ 1 M Q the halo mass function of our WDM three 
runs agree well with each other. The suppression with respect 
to CDM reaches a factor or two at M = 4 x 10 12 /i _1 M©, but 
it can be as large as a factor of 20. These values are consistent 
with prev ious studies, though sli ghtly stronger than those re- 
ported bv lSchneider et al.l (|2012h , however, the discrepancy is 
most likely due to the different halo definitions. The qualita- 
tive agreement with other works also is not surprising, given 
that the Tree-PM method has been the choice of those studies 
and it agrees with our other two methods. It is when we con- 
sider smaller masses, were this is no longer true, that we can 
enter into a regime hidden to previous simulations. 

Below M ~ 3 x 10 11 /i _1 M0, and for over an order of 
magnitude in halo mass, the T4PM and PM methods deviate sys- 
tematically from the Tree-PM run. The characteristic upturn 
in the halo mass function produced by artificial fragmentation 
is not present in either the PM or the T4PM run. The differ- 
ences reach a maximum factor of ~ 10 and 7, respectively at 
our mass resolution limit M~2x 10 10 /i _1 M©. All of this is 
consistent with the qualitative picture provided in Fig. [T] 

On the other hand, despite the lack of spurious fragmen- 
tation, there is no sign of a sharp cut-off, even in the T4PM 
run, as expectations raised from previous works suggest (e.g. 
iBenson et al.l 120131 ; ISchneider et al.ll2013l ). The abundance of 
low-mass objects only decreases slowly and shows a mild up- 
turn at rsj 5 x 10 10 H~ 1 Mq. We explore this issue in more detail 
next. 
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Figure 4. Images of randomly-chosen objects in six disjoint mass bins. These mass bins are equally spaced in log M over the mass range 
2 x 10 10 < M2oo/(h~ 1 Mo) < 5 x 10 12 . Each image displays the logarithmic projected density field computed using the T4PM method. The 
extent of each image is equal to 2 x R200, i.e. twice the virial radius of the respective halo. 



3.2 The nature of collapsed structures 

In order to explore the nature of the objects below the mass 
cut-off, and the origin of the mild upturn at ~ 5 x 10 10 H~ 1 Mq 1 
we have visually inspected all haloes found above our resolu- 
tion limit in the T4PM run. In Fig. U we provide density projec- 
tions of six randomly chosen objects found in six disjoint mass 
bins, which serve as examples of the type of objects that pop- 
ulate different regions of the halo mass function. We display 
mass carriers in a region of twice the size of the virial radius 
around each target halo. 

It is readily apparent that truly different structures are 
identified at the different mass scales. We now enumerate the 
most common features found in different mass bins: 

1) In the leftmost column (smallest mass bin), we find den- 
sity peaks just undergoing collapse, with highly disturbed mor- 
phologies, irregular boundaries and that usually show no clear 
center. In addition, we also find locally overdense regions that 
are however incompatible with the concept of a virialized dark 
matter halo. Most of these correspond to caustics - usually lo- 



cated at the radius where the particle orbits first turn back 
after crossing a potential minimum. Another, somewhat less 
common, occurrence of these are the centres of very dense 
filaments, and also the caustics of filaments. 

2) In the next mass bin, we typically find objects where 
there is clearly one of the three axes that has collapsed re- 
cently. These objects usually show a roughly round exter- 
nal iso-density contour, and a bar-like feature at their centre, 
which is the remnant of the filament whose folding produced 
the collapse of the objects. 

3) Objects in the third mass bin show less strong distur- 
bances. They correspond to roughly spherically symmetric ob- 
jects, but they clearly show many caustics, resulting from the 
continuous folding of the phase-space sheet. Commonly, they 
also show a bar, as those in the previous mass bin. 

4-5-6) Finally, in the three most massive bins, we find systems 
similar to those we usually find in CDM simulations and that 
can be unequivocally categorised as fully collapsed DM struc- 
tures, with a well defined centre and approximate spherical 



© 2013 RAS, MNRAS Q00.fTHT2] 



8 Angulo et al. 



10 1 



10 1 



© 



10 1 



10 1 



10 1 



: ' ' ' 


■ ' "i 




iU'" | ' ' 'd" 


1 ■ CO ■ ' 


1 IMJ 








oo oO o 


° O o 


/ : 


o 


o o 


D 


S °» 'oo 

O CO 


°y 


x : 


■o6> 
° o 

_ o o 


<b 




o %£ OOc 




■ 


6> ° 




o° o^%pJk 


















V 

o o o 


8 ^°° 

o Vl 


o 

1 


O 


■ 


- ° o 








- 


l°Vo 


c?i?> 








; 




^yf^ J^' 








■ 


C/* *■ A^ 








■ 


DO (A C2>-| 


^M>'P 7 









■ 


^fe«5 


cP o/^ 










bo o u oCO ^ 








: 


o o 












X) o 


/uo 










3 ° r, / 












■0 feA 


o 

o 








■ 




o o 








. 


% . •. . 


....I 


i 


1 1 1 


1 . I I . 


1 1 1 1 



10 1 



10 1 



10 1 



10 1 



10 1 



(z=0.5) [h" 1 M.] 



formation O However, haloes of different initial masses grow by 
significantly different amounts. Those of 10 13 H~ 1 Mq increase 
their mass by 30%, on average. At the expected half- mode 
mass scale, 2 x 1O 12 /i -1 M0, the increment is typically a factor 
of 2. Whereas, at 10 11 h~ 1 M Q it is a factor of 15! This very 
rapid mass increase at low masses has the consequence that 
most of the haloes below 10 11 h~ 1 MQ (e.g. those in the first 
two columns of Fig. |4j) are well above this limit by z = 0. This 
is consistent with the picture given above in which the objects 
found below the mass scale are simply a transient stage of halo 
formation, thus they quickly increase their mass and sit above 
the cut-off scale, where the growth proceeds at a slower pace. 
This is to be contrasted with the CDM case, where low-mass 
haloes have the earliest formation redshifts. 

Out of the 1413 points we display, there are 47 which 
are located below 1 x 10 11 /i _1 M© at both z = 0.5 and z = 
0.0. These haloes are not compatible with our interpretation, 
and they could be rare occasion where our mass resolution is 
not sufficient to avoid absolutely all fragmentation. Despite 
this, this is a very small population, which will not affect our 
results. 



Figure 5. The relation between halo mass at z = 0.5 and its 
descendant at z = 0. The red solid line indicate the median descen- 
dant halo mass at z = for ten bins of width Alog 10 M = 0.5 in 
progenitor mass at z = 0.5. Open circles show results for all the 
haloes at z = 0.5, whereas orange symbols highlight only the most 
massive progenitor of z = halos. The diagonal black line denotes 
a 1:1 relation. 



symmetry. The objects have much more clearly undergone an 
isotropic virialization than objects on lower mass scales. 

It is interesting to note that the evolutionary state of 
structures of increasing mass resembles different stages of halo 
formation in Hot dark matter cosmologies. Firstly, the lo- 
cal tidal field halts the expansion of a density perturbation 
along one axis, which eventually collapses, creating a 'pan- 
cake'. Then, collapse along a second orthogonal axes produces 
a filament. Material is accreted along that filament and the 
third and final axis collapses. Then, a relaxation process oc- 
curs, which finally gives rise to a DM halo in the sense of an 
approximately spherical and virialized object. Concurrently, 
the winding of the phase-space sheet continuously creates new 
caustics. Initially, these caustics are strong features, but as 
new ones are created and relaxation processes take place, they 
blend to produce a smooth density profile. If the connection 
between formation stage and halo mass is true, our picture 
suggests that haloes below the cut-off scale are simply haloes 
during early stages of their formation. This connection is also 
rather natural, as haloes at the cut-off scale should be forming 
monolithically as in HDM, not from accreting smaller objects 
as in CDM. 

This picture is supported by Fig. [5] where we show the 
descendant halo mass at z — of haloes detected at z — 0.5. 
We link haloes at these two redshifts by finding the object at 
z = that contains the majority of the particles associated 
with a FoF halo at z — 0.5. We highlight as filled orange sym- 
bols the z = 0.5 haloes that are the most massive progenitor 
of a z = halo. 

From this figure we can see that almost all haloes increase 
their mass consistent with a hierarchical picture of structure 



3.3 The abundance of virialized structures 

With the ideas discussed above in hand, we now return to the 
issue of the halo mass function in WDM cosmologies. Upon 
visual inspection of the members of our FoF- SO catalogue, it 
became obvious that many of the entries did not comply with 
the features usually found in halo catalogues built from CDM 
simulations. Consequently, we visually inspected and classified 
all haloes in our T4PM run into one of three groups. 

1. "Not Halos": In the first category we include all objects 
that appear as clear failures of our halo finder algorithm. These 
enclose mainly three cases: one corresponds to outer caustics 
of large haloes (which sometimes are located further than the 
virial radius) , another to descendants of haloes that have flown 
trough a more massive system (these are stripped of most of 
their mass, but their core survives). The third case correspond 
to dense sheets and filaments where sometimes the collapse of 
a further axis has started. 

2. "Proto-Halos" Our second category contains haloes 
that are not fully formed yet, but show clear isolated 3D den- 
sity enhancements. Here all three axes have collapsed, but the 
density peak has not fully virialised: we include here all ob- 
jects from highly anisotropic systems, that appear just after a 
violent collapse, to much more quiet haloes, where only minor 
departures from a smooth mass distribution exist. 

3. "Halos" The third and final category contains systems 
which can be unambiguously defined as a halo in the tradi- 
tional sense of approximately spherical objects showing clear 
three dimensional virialization, and that resemble those seen 
in CDM simulations. 

We note that we attempted to perform automatic classifi- 
cations using several different halo properties. Unfortunately, 



4 We find that the small number of haloes that see a reduction of 
their mass are systems that were accreted by a larger halo, their 
outskirts removed by tidal stripping, but the denser core survives 
in a orbit that yield them to outside the virial radius of the host 
halo, and thus are identified as separate haloes but with a reduced 
mass. 
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Figure 6. Relative contribution of different types of objects to 
the total WDM halo mass function. The black histogram shows 
the fraction of "Halos" . The green histogram shows the abundance 
of "'Proto-Haloes', whereas the green line indicates the fraction of 
objects that our SO-FOF algorithm wrongly identified as haloes. 
See the text for more details about our classification method. 



none of them could satisfactory separate the cate gories men- 
tioned above. Some of the measures introduced in lAbel et al.l 
([2012n or methods inspecting the shape of the Lagrangian 
pat ch of the identifie d structures (such as that proposed by 
e.g. iLovell et al.l 120120 may help in automating such a proce- 
dure. On the other hand, since proto-haloes, as well as some 
"Not Halo" objects, are likely simply early stages of halo for- 
mation, they are also likely to correspond to peaks in the initial 
conditions with only their collapse time differing from those 
corresponding to "Halos" . Thus, Lagrangian approaches might 
not clearly separate our three classes of objects. An additional 
complication for automatic classifications is that the haloes in 
the critical regime, < 10 11 h~ 1 MQ 1 are resolved with only a 
few thousand particles which is not enough to perform a de- 
tailed analysis of their internal structure. We will defer further 
exploration of these issues to future work. 

In Fig. [(J we show the relative contribution of each of 
these three categories to the WDM halo mass function. It is 
very interesting to see that the groups are clearly localized 
at different mass scales, although some overlap exists. High 
masses are dominated by standard haloes. Right below the 
cut-off, recently collapsed systems dominate. And the lowest 
masses receive a similar contribution from Proto-Haloes and 
from failures of our FoF-SO algorithm. 

Before we continue, we would like to note that, as in most 
classifications, the division between these three groups is some- 
what arbitrary. This is accentuated by the subjective nature of 
our visual inspection. For these reasons, we emphasise that the 
distinction between different categories just provides a qual- 
itative assessment of the nature of objects at different mass 
scales, and of how they affect the WDM halo mass function. 

Another point to note is that the fine division between 
the categories does depend on the force resolution employed. 
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Figure 7. Contribution of different types of objects to the WDM 
halo mass function. The red line show the abundance of standard 
dark matter haloes. The green line represents haloes in final stages 
of formation, while the blue line displays the abundance of objects 
in initial stages of formation. Finally, the magenta line shows ob- 
jects incorrectly identified as haloes by our algorithm. See the text 
for more details on the classification, and Fig. 0] for examples of 
structures in the various categories. We also display the mass func- 
tion expected for a CDM case for comparison (black line). Vertical 
dashed line indicate a limit where the abundance of haloes is not 
affected by finite force resolution. The bottom panel shows the ratio 
of our results to the expectations in a CDM scenario. 



We have explored this by carrying out our T4PM run using a 
PM mesh a factor of eight smaller, thus degrading our spatial 
resolution by a factor of two. There are three aspects worth 
noting. 

1) The amount of filaments/sheets in our catalogues, as 
well as the sum of Haloes plus Proto-Haloes remains roughly 
the same when the force resolution is varied. This is because 
the time of collapse of a filament depends mostly on large-scale 
tidal and density fields, which are less sensitive to the force 
resolution. 

2) The distinction between Haloes and Proto-Haloes is 
very different at different force resolutions. With higher force 
resolution, caustics are created more rapidly, there is more 
mixing, and haloes appear more relaxed. Note that due to 
computational limitations, it is not possible for us to increase 
the mass resolution of our runs needed to increase further the 
force resolution. 

3) The frequency of some type of FoF failures changes 
considerably with force resolution. In this case, we find that 
the number of Not-Haloes at the low mass end increased sub- 
stantially, mainly due to an increase in the number of caustics 
- a lower force resolution allows the turn-around radius to 
move outwards. 
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From this, we can conclude that the mass functions for 
the Haloes and the Haloes plus Proto-Haloes samples should 
provide the range in which we expect the mass function of vit- 
alised haloes to lie. 

We are now in the position to provide the most important 
result of our paper. In Fig. [71 we show the WDM halo mass 
function for two catalogues. The first one, denoted by a green 
dotted line, shows the abundance of systems categorised as 
Haloes, the second case, denoted by a blue dashed line, adds 
in the contribution of Proto-Haloes. Thus, the blue and green 
lines provide upper and lower limits to the abundance of col- 
lapsed objects in WDM cosmologies. These data are well fitted 
by the following functional form: 



^cdn 



-(M) 



1 + 



Mi 



1 + erf ( log 



MX 

~M~ 2 ) 



(6) 



with M\ — 3.9 x 1O 12 /i -1 M0 set by the half-mode mass scale, 
and M2 corresponds to the location of the small-scale cut- 
off, which we find to be 5.2 x 10 10 for the Haloes catalogue, 
and and 2.1 x 10 11 h~ 1 M Q for the Haloes plus Proto-Haloes 
sample. The best fits are displayed as dashed lines in Fig. [71 

We clearly see that once we neglect the contribution of 
failures of our halo finder, the WDM halo mass function shows 
a strong cut-off at small masses, and the upturn seen in Fig. [3] 
essentially disappears. The strong cut-off implies that there 
are no collapsed peaks below some scale, in agreement with 
what one would naively expect from the cut-off in the transfer 
function which predicts no small-scale density perturbations. 
This rules out the scenario in which a substantial population 
of low-mass haloes are created as a result of nonlinear self- 
interactions of the density field. However, a quantitative com- 
parison between the cut-off in the transfer function and that 
in the mass function is not straightforward, as the substantial 
differences between our two samples indicate. 

To illustrate this point, in Fig. [71 we also show two predic- 
tions for the shape of the WDM halo mass func t ion. D otted 
line shows the recent results of ISchneider et al.l (|2013r ) that 
is based on a sharp-k filter calibrated with a Tree-PM simula- 
tion, w hereas the solid line is the prediction of lSchneider et al.l 
(J2012|) based on the EPS formalism and extrapolating the re- 
sults of N-body simulations. Both predictions differ largely. 
And although the former prediction roughly matches the 
Haloes sample, as we discussed below, this sample only pro- 
vides a lower limit for the halo mass function. The later pre- 
diction, on the other hand, does not account for any cut-off in 
the halo mass function. 

The lack of an unique answer is a natural consequence of 
the complexity of structure formation and of the need for a 
common definition of a "halo" when numerical simulations and 
analytical formalism are compared. This is also true in CDM, 
but accentuated in WDM since low-mass scales are dominated 
by the collapse of filaments and by systems in the process of re- 
laxation or collapse. Despite the added complexity, examining 
the cut-off scale can be extremely useful to isolate successful as 
well as problematic aspects of analytic formulations, which in 
turn could lead to improvements in their foundations and ul- 
timately to a better understanding of structure formation and 
assembly in all cosmological models of structure formation. 



4 CONCLUSIONS 

In this paper, we were able to overcome the notorious difficul- 
ties associated with numerical simulations of the evolution of 
warm dark matter models. For decades, the numerical frag- 
mentation of filaments created an artificial population of low 
mass haloes which dominated the halo mass function on small 
mass scales. Here, we demonstrated that it is possible to avoid 
such artefacts by employing a force resolution consistent with 
the mass resolution of the simulations. With this, and for the 
first time, we could explore the halo mass function below the 
cut-off scale. We discovered a picture more complex than that 
present in CDM simulations. 

Structure formation in scenarios that have a small scale 
cut off in the power spectrum proceeds quite differently than in 
the CDM model. In CDM, haloes form mainly from accreting 
smaller haloes, which have the earliest formation times and 
the longest time to relax. On the contrary, in WDM, a large 
fraction of haloes forms from the direct collapse of filaments 
and small-mass haloes are typically the most recently formed, 
and thus are still in a process of virialisation. 

We find that this formation mechanism is imprinted in 
the WDM halo mass function: essentially different stages of 
halo formation dominate the counts at different masses. On 
the smallest mass scales we can resolve, identified objects are 
typically centres of filaments that are starting to collapse. On 
intermediate mass scales, objects typically correspond to fluc- 
tuations that have collapsed and are in the process of relax- 
ation, whereas the high mass end is dominated by objects 
similar to haloes identified in CDM simulations. 

In addition, we found that traditional group-finders pro- 
duce, on small mass scales, catalogues with objects not consis- 
tent with the definition of a halo. In a CDM calculation, essen- 
tially all dense structures are always part of haloes, and even 
the simplest approaches, such as the FoF algorithm, succeed in 
selecting appropria t e objects above a certain mass limit (e.g. 
IWarren et al1l2QQ6l : iMore et al.ll201lh . On the other hand, in 
WDM these approaches prove to be unreliable. The small-scale 
cut-off in the primordial transfer function, implies that there 
are dense structures (e.g. filaments, caustics of large haloes) 
that can be incorrectly considered as haloes. These misidenti- 
fications dominate the WDM mass function below the cut-off 
scale. For the simulations we considered here, we could by- 
pass this problem by visually inspecting and classifying differ- 
ent systems. This is not prohibitively demanding, thanks to 
the relatively low number of systems in our runs. However, of 
course, this is in general not true, and improved halo finders 
are desirable for future work. 

After neglecting these failures of halo identification, we 
observe a strong cut-off in the halo mass function, with very 
few objects found below the cut-off scale. These correspond to 
haloes undergoing rapid collapse and virialisation. Our results 
indicate that the cutoff scale in the initial density fluctuations 
does indeed translate into a comparably strong cutoff scale 
in the halo mass function (Fig. [6}. This implies that filaments 
and sheets that formed in early stages of gravitational collapse 
are remarkably stable and do not fragment due to the lack of 
small scale perturbations. 

Our work poses several questions. The most natural one 
is the role of artificial fragmentation in the internal properties 
of DM haloes. In the absence of spurious haloes, filaments 
get denser and the accretion of material happens continuously 
from filaments, not as a sequential accretion of small dense 
haloes which would likely behave differently dynamically. It is 
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not clear whether this difference will impact, for instance, the 
concentration of haloes near the cut-off scale. Thus, we plan 
to carry out a suite of higher-resolution simulations to address 
this. 

Another question arises from the fact that the exact mass 
scale at which the cut-off in the halo mass function is located, 
does depend on the exact definition of what constitutes a 
halo. Especially, it depends on the separation between a fully 
formed halo and a halo in the process of formation and/or 
relaxation. A priori it is thus not clear what one may wish 
to call a halo. Functional definitions might involve regions of 
sufficiently high DM density that would allow the associated 
baryons to collapse and cool, or surface densities that allow for 
gravitational lensing or more theoretically inspired measures 
such as a limit on the velocity anisotropy of the dark matter 
particles or perhaps particular axis ratios inferred from the 
velocity ellipsoid or inert ial tensor. Whatever the exact defini- 
tion one may choose the abundance of "haloes" fulfilling these 
choices will likely vary quite dramatically. 

Therefore, any key ingredient for a new halo finder will be 
related directly to the properties that one aims to select for. 
For instance, the isotropy of the velocity dispersion, or the 
number of streams, an d include the requ irement that all three 
axes have collapsed (|Falck et al.l l2012h perhaps augmented 
with l ocal measures of t he velocity dispersion in all three direc- 
tions (|Abel et al.ll2012r ). A common definition is also needed 
for a proper comparison with analytic models for structure 
formation. It will be desirable for future work to investigate 
whether there are generic definitions that sensibly define a 
concept of dark matter halo that is applicable to CDM as well 
as WDM computational cosmology questions. 

Even once the issue of halo definition is decided upon, 
there is another problem related to the numerical methods. 
We find that the exact virialisation state as well as the number 
of caustics depend sensitively on the force resolution employed 
in our simulations. This is different from standard CDM runs 
because of the difference in the formation mechanism, and 
because of a special combination of force and mass resolution. 
However, this also opens an exciting possibility: because of 
the free-streaming scale in WDM, there is also a upper limit 
to the density and also a scale for small-scale gravitational 
interactions, which could eventually allow us to simulate the 
full range of length and mass scales relevant for the formation 
of a halo. 

The question of halo definition, and the role of artificial 
fragmentation, is closely related to the issue of what is the 
minimum halo mass and/or virialisation stage for a halo to 
host galaxies. Fragmentation of gas into stars may well still 
occur in places not regarded as haloes in the traditional sense, 
perhaps already in regions of filamentary collapse: anywhere 
where a local gravitational potential well is already aggregat- 
ing the baryons. On the contrary, the continuous folding of 
the phase-space sheet might continuously shock-heat the gas, 
and the cut-off in the respective galaxy stellar-mass function 
might not trivially relate to that in the halo mass function. 

Unfortunately, given the strong limits on the mass of 
the potent i al WDM particle from Lyman- a forest constraints 
IViel et al.1 d2Q05h : iBovarskv et~ai1 (J2009aT ) the mass scale on 
which such a different galaxy formation scenario could be re- 
lated to the observable Universe is severely limited, and other 
physics might be more relevant. Nevertheless, WDM scenar- 
ios remain of great interest to obtain a much better theoretical 
understanding of how dark matter haloes assemble and how 



the collisionless fluid virializes. In particular, they very closely 
resemble the events that lead to the very first dark matter 
haloes even in a CDM s cenario albeit on rad i cally different 
mass and length scales ( Diemand et al.l 120051 ; ICoerdt et al.l 
120071 ; llshivama et alJuQlOl ) . We foresee our results stimulating 
a more detailed exploration of the formation of the smallest 
structures expected to form in a given cosmology, which will, 
hopefully, advance our overall understanding of structure for- 
mation. 
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